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ABSTRACT 

Least-squares frequency switching (LSFS) is a new method to reconstruct signal 
and gain function (known as bandpass or basehne) from spectral line observations using 
the frequency switching method. LSFS utilizes not only two but a set of three or 
more local oscillator (LO) frequencies. The reconstruction is based on a least squares 
fitting scheme. Here we present a detailed investigation on the stability of the LSFS 
method in a statistical sense and test the robustness against radio frequency interference 
(RFI), receiver gain instabilities and continuum sources. It turns out, that the LSFS 
method is indeed a very powerful method and is robust against most of these problems. 
Nevertheless, LSFS fails in presence of RFI signals or strong line emission. We present 
solutions to overcome these limitations using a flagging mechanism or remapping of 
measured signals, respectively. 

Subject headings: methods: data analysis — techniques: spectroscopic 



1. Introduction 

Radio astronomers were quite fast from the 
beginning confronted with problems as un- 
known intermediate frequency (IF) gain func- 
tions and instabilities, both in frequency and 
time, which makes spectral analysis of faint 
or very broad emission line sources very com- 
plicated. Up to now mainly two different 
schemes are used to overcome these problems: 
position switching and frequency switching. 
Both methods reduce the receiver gain insta- 
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bilities by measuring a reference spectrum, ei- 
ther off-source (position switching) or with 
a detuned local oscillator (LO) frequency to 
move the line of interest outside the observed 
spectral band (frequency switching). There 
are a lot of drawbacks applying both meth- 
ods. Using position switching it is necessary 
to avoid different continuum levels (or very 
broad line emission) towards the on and off 
position which of course limits the usage of 
position switching in Milky Way observations, 
having emission from all directions. This is es- 
pecially valid for the important spectral lines 
of neutral atomic Hydrogen or CO, due to 
their high area filling factor. Furthermore, 
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there is loss of observing time either when 
redirecting the telescope between on and off 
positions or during the retuning of the LO 
frequency, respectively. Frequency switching 
suffers also from the gain curve changes while 
shifting the spectral range. The most prob- 
lematic aspect of both methods is the loss of 
half of the observing time. In-band frequency 
switching (both LO phases provide the line 
of interest lying in the observed bandwidth) 
avoids this loss. Unfortunately, one loses half 
of the available bandwidth (which is equiva- 
lent to the loss of velocity coverage or number 
of spectral channels), because a proper sep- 
aration of both signals is needed. Observ- 
ing single objects this is usually not a ma- 
jor problem, but when performing, for exam- 
ple, a blind survey one would strongly suffer 
from such a restriction. We also point out, 
that many important questions of modern cos- 
mology rely heavily on highest-sensitivity ra- 
dio measurements of the high-redshift uni- 
verse. Modern spectrometers based on field 
program mable gate arrays (FPGA) technol - 
ogy fe.g. lBenz etaPboOSl : IStanko et al.ll2005l ) 
are best suited for that purpose but meaning- 
less without sophisticated methods for band- 
pass calibration. 

Heiled ^200i ) presented a new methoc. 
called Least-Squares Frequency Switching 
(LSFS). LSFS is able to deal with all problems 
discussed above, making it the best choice 
for future spectral line observations with ra- 
dio telescopes especially well suited for H L 
However, it requires minor hardware changes 
at the telescope in order to provide not only 
two, but a set of three or more LO frequencies 
within one switching cycle. This is not a sub- 
stantial problem, and there already is a work- 
ing sy stem using LSFS at the Arecibo tele - 
scope (|Heilesll2007l : IStanimirovic et al.ll2006l ). 

Of major interest in radio spectral observa- 
tions — especially in high-redshift H I astron- 
omy — is the ability of a "bandpass removal" 



tool to provide high-quality results in a statis- 
tical sense. Most observations today need to 
integrate at least a couple of minutes to reach 
the desired sensitivity limit. While in theory 
the noise level scales as 1/ \/ tAu according to 
the radiometer equation, with Av being the 
bandwidth and t the integration time, this is 
not necessarily true for a real receiving sys- 
tem. Modern backends have Allen times, to, 
(for t < to the radiom eter equation holds ) 
of hundreds o f seco nds ( Stanko et al.l 20051 ). 
Winkel et all tooil ) show that in presence of 



radio frequency interferences (RFI) this can 
be limited to less than few tens of seconds. We 
analyze, if the LSFS method has an impact on 
the RMS level (the sensitivity) achievable. It 
turns out, that in statistical sense the LSFS 
performs very well, yielding only a slightly de- 
crease of sensitivity. We also test the robust- 
ness of the LSFS versus several typical prob- 
lems at radio telescope sites. These are for ex- 
ample RFI events, possibly bandpass instabil- 
ities in time and frequency, as well as (strong) 
continuum sources. LSFS works well under 
most tested circumstances except for RFI and 
strong emission lines. However, small changes 
to the original LSFS method already provide 
meaningful workaround mechanisms. These 
are discussed in subsequent sections. 

Section [2] contains a br ief suinmary about 
the LSFS method based on'Heilej (|2007l ). The 
robustness of LSFS is analyzed in Section [3l 
In Section [4] we finally calculate the comput- 
ing time needed in order to perform the vari- 
ous calculations related to the LSFS method. 
We propose to use a modified method for the 
computation of the SVD, based on an algo- 
rithm available for sparse matrices. This im- 
proves computation speed for the (one-time 
per LO setup) calculation of the correlation 
matrix by an order of magnitude. Section [5] 
contains our summary. 
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2. Least-squares frequency switching 
(LSFS) 

Most spectroscopic observations in radio 
astronomy use the heterodyne principle where 
the radio frequency (RF) signal is multiplied 
with a monochromatic signal of a LO. An ap- 
propriate low-pass filter applied after this op- 
eration provides the desired IF signal at much 
lower carrier frequencies. The whole system 
can be described by 



^'if(/if) = G'if(/if)Grf(/rf)x 
[TA{fRF) + TA + TR{fnF)+TR] 



(1) 



with PiF being the power in the IF chain. Gjf 
and Grf sue the gain functions at IF and RF 
stage, respectively. The gain acts on the sig- 
nals which enter the feed — the astronomical 
signal of interest plus the contribution from 
the sky which we denote as Ta — as well a s 
on the noise of the receiver, Tr. HeilesI ( 2007 ) 
separates Ta and Tr into a frequency depen- 
dent and independent (continuum) part. Note 
that the gain is not simply a scalar but has a 
spectrum due to the filter curves. 

To recover from the measured signal, Pip, 
the signal of interest, Ta, one needs to know 
the gain spectrum Gif (Grf can be treated as 
constant with frequency on modern receivers). 
Traditionally, this is achieved by measuring a 
reference spectrum without any spectral fea- 
tures either by position or fr equency switch- 
ing. We refer to iHeilesI (120071 ) for a review of 
position and frequency switching. 

Equation ([T]) reads in a simplified form as 



-Pif(/if) = Gif(/if)'S'rf(/rf) 



(2) 



where we combined all signals enteri ng the 
mixer to 5'rf(/rf)- In contrast to iHeiled 
(|2007l ) we dropped the assumption that the 
input signals superposition of frequency 
dependent and independent parts. Using 
modern broadband spectrometer backends 
we likely can not treat the continuum signal 



as constant over the entire observed band- 
width. Moreover we might also be interested 
in the continuum emission itself. Broadband 
spectrometers could also be used to gener- 
ate continuum maps. Furthermore, from the 
mathematical point of view it is not necessary 
to perform the separation — all subsequent 
equations do not rely on it. 

As before, we are interested in obtaining 
Gif(/if)- Introducing not only 2 but a set 
of different LO frequencies (for a detailed 
analysis of h ow to choose appropriate LO fre- 
quencies, see 



HeilesI 120071 ) . we end up having 



A^ • I equations 



i+Ain ■ 



(3) 



In this representation we now use integer in- 
dices representing the spectral channels of the 
backend. Then i is the i^^ out of / channels, 
Ain is the frequency shift of LO n versus LO 
(the unshifted LO) given in channels. By us- 
ing different LO frequencies we of course ob- 
serve somewhat different spectral portions of 
the input spectrum S. Without loss of gen- 
erality we normalize the input signal to have 
a mean value of unity, Si+Ain = 1 + Si+Ain > 
leading to 



Pi,Ai„ — Gi + GiSi^Ai 



(4) 



which can be solved using nonlinear least- 
squares techniques. However, iHeilesI (j2007l ) 
converted the equation to an iterative linear 
least-squares problem by solving for the differ- 
ence of guessed values of Gf and sf+Ai„ ^o™ 
their true values. From these guessed val- 
ues one can of course compute the associated 
output P^Ain power for each spectral channel 
and LO setting. After some si mplifications 
(dropping higher order terms; see lHeilesll2007l ) 
equation dH) transforms into 
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i,Ain 



Gf 



5Gi_ 



+ 6s 



i+Ai„ 



(5) 
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The (5-terms denote the difference between the 
true and the guessed value of the correspond- 
ing quantity. A further constraint is needed 
in order to keep the mean RF power approxi- 
mately constant, namely X^i n "^^i+Ain = 0. 

For convenience we use matrix notation for 
equation ([5]) 



P = 
= (pIo, 



Xa 



(6) 
(7) 



a 



= [go, 



T 

'Pi.N-1 



,<5s/-l+Aijv_i)(8) 
9i = — • (9) 



Least-squares fitting is achieved by computing 
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a . 



a 



(10) 



with the covariance matrix a. Computing 
a requires matrix inversion which in general 
does no t exist necessa rily. To deal with degen- 
eracies, Heiles (|2007l ) proposes the Singular- 
Value Decomposition (SVD) of matrix X 

X = U[W]V^ (11) 

with the diagonal matrix W containing the 
so-called singular values Wi. In the case of 
degeneracies one or more of the Wi are close 
to zero, leading to infinite (or huge) numbers 
when inverting. It turns out that 



(aXT) = V 



1 

W 



(12) 



The critical singular values can be treated 
separately (e.g. setting the inverse values to 
zero). By computing the SVD of the matrix 
X one can directly solve equation ^ without 
encountering any problems caused by degen- 
eracies. The computation of the SVD of a ma- 
trix, e.g. for N = 8, I = 1024, is possible on a 
modern PC but is not finished within fractions 
of a second; see Sect. d] for details. Neverthe- 
less, the SVD calculation fortunately needs to 



be done only once per LO setup, as the matrix 
itself is independent from the measurements. 

The assumed normalization of the signal 
seems to be a somewhat arbitrary assumption. 
But if we assume the bandpass (gain) curve 
to be normalized, we can in practice easily 
attribute associated gain factors to the sig- 
nal which nevertheless has to be calibrated in 
terms of intensity. This way we can uniquely 
reconstruct the overall power of the input 
signal (in arbitrary units) by computing the 
mean, m, of the measured signal, normal- 
izing (dividing by m), calculating the LSFS 
which gives a signal of mean value of unity 
and finally multiplying the reconstructed sig- 
nal with m. It is clear that this scheme will 
only work if the gain curve remains constant. 
This can be expected at least for the duration 
of the observing session so that the computed 
gain factor remains constant for a single ob- 
servation. 

The LSFS method also allows reconstruct- 
ing the continuum part of the input signal. 
While for position and frequency switch- 
ing the separation of RF-dependent and - 
independent parts was necessary, the LSFS 
algorithm allows the reconstruction of the 
complete mixture of signals which are put into 
the mixer. Of course this implies that there 
is possibly the need for further disentangling 
these signals into line and continuum com- 
ponents (from astronomical sources, ground, 
and receiver noise) which may even have dif- 
ferent spectral indices. 

The equations hold for small values of 
Si+Ai„ as they were computed only to first 
order approximation. In most cases this is 
easily fulfilled in radio astronomy as the ob- 
served lines are much weaker than the typical 
intensity of the unavoidable continuum level 
produced by the atmosphere, ground, and 
receiver noise, which sum is known as the 
system temperature Tgyg. However, there is 
one case known where we indeed have a sig- 
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nal much brighter than the continuum level: 
the H I emission of the Milky Way which 
can reach intensities of a few 100 K while the 
system temperature for a typical telescope is 
20 — 40 K. We will address this issue in the 
following analysis. 

3. Robustness of the LSFS 
3.1. Setup 

We implemented the LSFS algorithm within 
the programming language C and performed 
various tests to investigate statistical stability, 
response to possible variations of bandpass 
shapes, impact of RFI signals, and its ability 
to deal with (strong) continuum sources. 

For our testing purposes we generated 
spectra (1024 spectral channels) by simulat- 
ing several Gaussian-shaped (faint) emission 
lines of different intensities and widths on 
top of a constant signal (which shall resemble 
those continuum signals with spectral index 
of zero). After adding Gaussian noise, these 
emission lines are partly not anymore visible; 
compare for example the signal spectrum of 
Fig. [2] and Fig.O This "true" signal is then 
multiplied with a gain function 



G'if(/if) 



■ ^iF ■ '-'IF y^"^) 



Ggt = i [tanh (5/ + 5) 



/^wave 
<-^iF 



G 



poly 
IF 



tanh (5/ -5)] (14) 
l + 0.1cos(F7r/) (15) 



which we adopted from Heiles ( 200?! ) to al- 
low for a better comparison of our results 
to that work. / is the frequency which we 
transformed to the spectral (i) domai n using 
/ = 2 .1(i - 512)/1024. In contrast to iHeiled 
^20071 ) we chose a spectral portion where the 
gain curve does not gets too close to zero. 
This would break the assumption that Gi is 
of order unity and distorts the normalization 
scheme which we pres ented in the previous 
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section. iHeiled (j2007l ) did not encountered 



Fig. 1 . — Raw input spectra as would be mea- 
sured by the receiving system using the MRS 
scheme. For better visualization the spectra 
(grey solid lines) were stacked and a noise-free 
analogon (black solid lines) was overplotted. 
Each spectrum is the multiplication of the 
"true" input signal and the IF gain function. 
Due to the different LO frequencies within a 
LO cycle the signals of interest are folded to 
different spectral channels. 

that problem because he neglects continuum 
emission. Note, that this is no drawback, as 
one can easily choose in practice those por- 
tions of the spectra which fulfill Gj ~ 1. We 
varied the two parameters around A = 0.1 and 
F = 4 to change the bandpass in amplitude 
and shape for some of our tests. A small vari- 
ation of F around 4 already has a dramatic 
impact on overall shape of the gain curve. 

3.2. Statistical stability 

We started by examining the statistical sta- 
bility of the method. For this purpose we 
generated spectra for a set of 8 LO freq uen- 
cies using the MRS scheme (|Heilesll2007l ): see 
Fig.[TJ As shown in Fig. [21 the solution for a 
single set of spectra does not necessarily pro- 
vides the "true" signal and bandpass shape, 
but there are small systematic effects. We at- 
tribute these partly to the influence of sin- 
gle (strong) noise peaks to the overall solu- 
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Fig. 2. — Reconstructed and original signal 
(top) and bandpass shape (bottom) of a sin- 
gle spectrum from our simulations. The signal 
is a superposition of noise and three unrelated 
line signals — two of them are well below the 
noise level but are visible after integration of 
several spectra; see Fig.[3] (top). Note that 
the noise level of the reconstructed signal is 
about a factor \/8 smaller because 8 spectra 
(the different switching phases) result in a sin- 
gle reconstructed signal spectrum. 



tion. Remember also, that the set of equations 
which we use to solve the decoupling of signal 
and gain is in linear-order approximation. It 
is important to note, that the iteration was 
in any case not interrupted until the solu- 
tion had converged. The question is whether 
these systematics cancel out after integration 
of several spectra or whether they remain. In 
order to have a measure of t he "goodness" 
of the solution we adopt from iHeilesI (|2007l ) 



two quality indicators — the RMS level of the 
reconstructed signal (denoted as RMS) and 
the RMS of the residual gain curve (denoted 
as a). The latter quantity uses the residual 
which is the difference between the true gain 
curve (noise-free) and the reconstructed gain. 
The first quality indicator, the RMS, is cal- 
culated making use of all spectral channels 
except those containing the signals of inter- 
est. For the purpose of comparison with the 
noise level of the signals we also rescale a with 
the gain factor, which was used to scale the 
signals. It is also obvious, that the recon- 
structed signals ideally should have a factor 
of \/8 lower noise compared to the originally 
generated signals. This is because each recon- 
structed signal was calculated using eight "ob- 
served spectra" (one LO cycle). Furthermore, 
we analyze the behavior of the indicators as a 
function of integration 'time'. This is done by 
successively summing up adjacent spectra. In 
each step this reduces the number of spectra 
by a factor two. Therefore, we start with a to- 
tal number of generated 1024 true spectra (or 
8 • 1024 measured spectra, respectively) which 
is a power of two. The LSFS was calculated 
for each of the 1024 spectra, then the sum- 
mation of the reconstructed signal and gain 
curve was performed stepwise. Theoretically, 
the functional dependence of RMS vs. inte- 
gration time is given by the radiometer equa- 
tion 



Pit) 



-0.5 



(17) 



\/Ai/ • t 

Fig.[3] clearly shows, that despite the fact 
that individual spectra were not perfectly 
handled the integrated signal as well as the 
bandpass visually match well. However, Fig.d] 
reveals increased RMS values of about 30% 
of the reconstructed signal and about 35% 
higher noise for a. Calculating the RMS and 
a values with respect to a 3rd-order polyno- 
mial (fitted after integration) results in sig- 
nificantly lower noise values which are only 
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Fig. 3. — After integration of 1024 subsequent 
spectra the faint signals previously hidden by 
noise show up. Both signal (top) and band- 
pass shape (bottom) are well recovered, as the 
functional behavior of the quality indicators in 
Fig. m reveals. 

slightly increased compared to the theoretical 
expectation value by 8% (RMS) and 10% (a), 
respectively. Obviously the residual system- 
atics can be described by a low-order polyno- 
mial. We point out, that the RMS behavior 
of the signal is of much greater interest from 
the observers point of view. If the gain curve 
is sufficiently stable with time, one can also 
compute the systems gain dependence to high 
precision by using thousands of spectra. The 
results show that one is effectively not losing 
sensitivity using this metho d as this w as the 
case in earlier attempts (e.g. Liszt 199?! ) with 
increased noise levels of about 100%. 
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Fig. 4. — Functional dependence of the dif- 
ferent quality indicators vs. integration time 
(scans). The boxes mark the noise level 
(RMS) of the true signal, the circles represent 
the noise of the reconstructed signal. The tri- 
angles (up) mark the RMS values of the gain 
residual, a, calculated using the difference of 
the true and reconstructed gain curves. After 
subtracting a third-order polynomial both the 
RMS (crosses) and a (triangles, down) quan- 
tities are closer to the theoretical value. 

3.3. LSFS and strong line emission 

In the previous section we addressed the 
possible problem of strong emission lines (as 
would be the case in galactic H I observations) 
which could violate the assumption of small 
variations of the (normalized) signal around 
unity. Here, we use a strong emission line 
to test its influence on the LSFS. To make 
it short — the LSFS fails completely; see 
Fig.O The intensity of the strong line signal 
at spectral channel 600 is about 5 times higher 
than the baseline level (system temperature). 
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Fig. 5. — In the presence of a strong emission 
line the LSFS method fails entirely. After nor- 
malization, all spectral features must be close 
to unity, otherwise the linear-order approxi- 
mation is no longer valid. A solution to the 
problem is remapping of the input signal; see 
Fig.Eland Fig.[71 

This causes heavy distortions during the re- 
construction process (comparing the relative 
amplitudes of the input and reconstructed sig- 
nal). 

We identified a workaround to the prob- 
lem: by remapping the observed (normalized) 
signal, P, in terms of a nonlinear function 
one can treat strong signals into the realm of 
small variations around unity. In our exam- 
ple we tried the mapping function P \/P, 
with x = 4. It is not clear, though, that 
the reconstructed signal and bandpass can be 
transformed back by simply using the inverse 
P P^. But indeed it turns out, that this is 
possible; see Fig.[6]and Fig.[71 

3.4. Bandpass instabilities 

To further test the statistical stability, we 
now change the bandpass shape and ampli- 
tude with time. First, we only changed the 
shape slowly but using the same shape for 
each bandpass within a single switching cy- 
cle (8 adjacent spectra have the same shape). 
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Fig. 6. — Strong emission lines (top) can be 
handled by remapping the measured signal by 
a nonlinear function, e.g. P — > ^/P. This en- 
sures the LSFS method to be in the linear 
regime. The bottom panel shows a zoom-in 
for better visualization. The quality indica- 
tors are shown in Fig. [71 

This should resemble the situation at the tele- 
scope site as we can (hopefully) expect the 
bandpass shape to be independent of switch- 
ing frequency, keeping in mind that the fre- 
quency shifts are very small compared to the 
total bandwidth. We could not find any sig- 
nificant difference to the undisturbed case; see 
Fig.Eland Fig.O 

As the slowly changing bandpass was no 
challenge for the LSFS algorithm, we also 
changed the bandpass shape more rapidly, but 
in a manner that there are no systematic dif- 
ferences between the different LO phases. The 
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Fig. 7. — Strong emission lines (top) can be 
handled by remapping the measured signal by 
a nonlinear function, e.g. P — > v^- This 
ensures the LSFS method to be in the lin- 
ear regime. The quality indicators (see Fig. [4] 
for the explanation of the symbols) show that 
the remapping works correctly in a statistical 
sense. There is no significant increase of the 
RMS or a values compared to the undisturbed 
case; see Fig.Hl 

outcome of this is shown in Fig. [TO] and Fig. 111! 
The LSFS method could not reconstruct the 
signal and gain curve. The residual is smooth 
but can only be described by a high-order 
polynomial. In fact, by computing the RMS 
with respect to a third-order polynomial we 
end up with significantly increased noise val- 
ues and a. 

For completeness, we also changed the 
shape of the bandpass in a systematic way 
by multiplying each bandpass with a lin- 
ear function which slightly drops off towards 
higher frequency (negative slope). The slope 
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Fig. 8. — A slow changing bandpass shape 
(constant gain curve during one LO cycle) has 
no measurable influence on the LSFS method. 
The signals can be treated as well recovered 
(see Fig. ED. 

of this function was steeper with higher shift- 
ing frequencies. This should mimic one of 
our early test observations with a digital 
fast-fourier- transform (DFFT) spectrometer 



proto type (jStanko et al.l l2005l : IWinkel et al 
20071 ). where we were forced to use an LO 



frequency far off any specifications. We com- 
bined this systematic error with a slowly over- 
all change of the bandpasses; see Fig. [12] and 
Fig. [13] This time, the outcome was slightly 
better — calculating RMS and a with respect 
to a 3rd-order polynomial leads to acceptable 
results in case of the signals RMS. However, 
the value of a does not decrease significantly 
after the summation of about 100 spectra. 
At the end the noise is about a factor of 
four higher than expected. The RMS level 
of the reconstructed signal is not significant 
increased compared to the undisturbed case. 
Without subtracting a baseline the RMS and 
a values are even nearly independent on inte- 
gration time. 

We note that the latter two cases of very 
strong bandpass instabilities are far from any 
realistic scenario at modern radio telescopes. 
IF filter devices may have response to temper- 
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Fig. 9. — A slow changing bandpass shape 
(constant gain curve during one LO cycle) has 
no measurable influence on the LSFS method. 
The quality indicators (see Fig. [J] for the ex- 
planation of the symbols) show no significant 
increase of the RMS or a values compared to 
the undisturbed case; see Fig.|H 

ature and frequency variations but on a much 
smaller scale than we used to test the robust- 
ness of LSFS against those instabilities. In 
case of slowly varying gain curves the LSFS 
performs as good as without bandpass varia- 
tions. 

3.5. Continuum sources 

When mapping a region of the sky one 
often encounters the situation that contin- 
uum sources contribute significantly to the 
observed signal. Using the typical in-band 
frequency-switching algorithms, we have to 
assume that the spectra of these sources are 
sufficiently flat not showing any signiflcant 
difference in the two switching phases. LSFS 
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Fig. 10. — A fast changing bandpass shape 
(see text) pushes the LSFS method to its lim- 
its. The bandpass and signal (top and bottom 
panel) were not well reconstructed. The dif- 
ference between true and reconstructed gain 
curve as well as the baseline of the recon- 
structed signal can not be described by a low- 
order polynomial; see Fig. [11] for the quality 
indicators. Note, that due to the rescaling 
of the signal the uncertainties are much more 
visibly prominent in the signal domain than 
in the gain curves. 

is much less dependent on this assumption as 
we switch only by a small fraction of the to- 
tal bandwidth. But the greatest advantage 
of LSFS in this context is that the continuum 
signal will be part of the recovered signal spec- 
trum, when using our normalization scheme. 
As spectrometer bandwidths have grown up 
to hundreds of MHz or even GHz nowadays it 
has become possible to also map continuum 
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Fig. 11. — Quality indicators (see Fig. H] for 
the explanation of the symbols) for a fast 
changing bandpass shape (see text). The dif- 
ference between true and reconstructed gain 
curve as well as the baseline of the recon- 
structed signal can not be described by a low- 
order polynomial. Both the RMS and a val- 
ues are much higher than in the undisturbed 
case and their functional behavior is far from 
linear. 

sources 'for free' within a spectroscopic obser- 
vation. 

Fig. dH shows the result for the case that a 
continuum source is superposed to the spec- 
tral lines. Its intensity is described by 



Iu = A{- 

.^0 



(18) 



with spectral index a = —2 and amplitude 
A = 2 assuming vq = 1420 MHz and a fre- 
quency resolution of 50kHz {5v ~ lOkms"^) 
per spectral bin. Both the continuum signal 
as well as the spectral lines were nicely recov- 



10.40 
10.30 
10.20 
10.10 
10.00 
9.90 
9.80 
9.70 
1.40 
1.20 
1.00 
0.80 
0.60 
0.40 
0.20 
0.00 





True Bandpass 
Reconstr. bandpass 



200 



400 600 
Spectral channel 



800 



1000 



Fig. 12. — In case of a systematic change of 
the bandpass shape which is due to the shift 
frequency (see text) LSFS fails to reconstruct 
the signal (top) and bandpass (bottom); see 
Fig. [13] for the quality indicators. 

ered, as the quality indicators (Fig. fTS]) show 
no increase in the RMS or a values. 

3.6. Radio frequency interference 

One of the key properties of each data re- 
duction pipeline used in radio astronomy to- 
day is the capability to deal with radio fre- 
quency interferences (RFI). These artificial 
signals are in general variable on timescales 
down to /is. Therefore, one of the most inter- 
esting analyses in this section is the impact of 
such interferences on the LSFS method. For 
simplicity, we started by adding two narrow- 
band interferences whose amplitudes obey a 
power law. This is — at least at the 100- 
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Fig. 13. — Quality indicators for a systematic 
change of the bandpass shape which is due to 
the shift frequency (see text). The residual 
gain curve can — to some extent — be de- 
scribed by a low-order polynomial, but after 
integration of 100 spectra the a value (see 
Fig. S] for the explanation of the symbols) no 
longer decreases. 

m telescope at Effelsberg — one of the most 
common types of interference. 

As the LSFS algorithm assumes the signal 
to be stable, it was not surprising that the re- 
sult is practically useless. Due to the coupling 
of channels with different frequency shifting, 
we end up with a number of contaminated 
spectral channels which is higher than the ini- 
tial number of channels affected; see Fig. [TBI 
The only solution is actually to address the 
RFI problem before performing the LSFS. 



Winkel et al 



(|2007l ) presented an algo- 
rithm which detects interferences down to the 
^ 4(Trms level. Having detected interference 
peaks, one can flag these data in order to 
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Fig. 14. — The presence of a continuum source 
does not have negative influence on the out- 
come of the LSFS algorithm. The signal con- 
tains a continuum source of spectral index 
a = —2 (top). Both spectral and continuum 
emission are well recovered as the quality in- 
dicators (see Fig. [T5]l reveal. 

exclude them from the computation. Flag- 
ging data points is equivalent to projecting 
the correlation matrix in Eq. d?]) to a subspace 
which does not contain contaminated spectral 
channels. This, however, would require to re- 
compute the SVD of the matrix each time 
the spectral channels containing RFI would 
change. This is far from practical as comput- 
ing the SVD for 1024 spectral channels and 
8 LO frequencies takes at least a few minutes 
on a modern PC. 

By far easier is the following alternative: 
setting all spectral channels containing an RFI 
signal (those are of course different channel 
numbers for different shifting frequencies) in 
p to zero. Of importance is here a robust cal- 
culation of the mean signal strength by drop- 
ping all disturbed spectral channels. Other- 
wise, the gain factors would depend on the 
actual strength of the RFI signals. 

As our algorithm is not able to find RFI sig- 
nals hidden in the noise (though an iterative 
scheme may be possible, were one performs 



12 



Quality indicators with integration time 



0.1 
w 

0.01 



0.001 



True RIVIS m 

Reconstr. RIVIS o 

Reconstr. RMS, base3 + 

Sigma a 

Sigma, baseS v 



1 10 100 

Integration time (number of spectra) 



1000 



Fig. 15. — Quality indicators (see Fig. H] for 
the explanation of the symbols) are not sen- 
sitive against continuum sources in the data. 
There is no significant increase of the RMS or 
a compared to the undisturbed case. 

the search for interferences at different inte- 
gration levels), we only set spectral channels 
to zero which contain an interference signal of 
> 4(Trms- We added three narrow-band RFI 
signals whose amplitudes obey a power law 
with spectral index v = —1.5. The leftmost 
signal at spectral channel 250 was persistent 
in all LOs except 1 and 2. The signal at chan- 
nel 600 was only added in every second LO 
and the rightmost interference at channel 680 
was added for each LO. 

The outcome is shown in Fig. [18] and 
Fig-dHl The bandpass was well recovered. 
However, each RFI leaves behind some 'finger- 
print' in the reconstructed signal, the residual 
strength of which obviously depends on the 
number of affected LOs. These "left-overs" 
are nevertheless easy to handle as the spectral 
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Fig. 16. — RFI signals can have a severe ef- 
fect on the solution of the LSFS. We added 
in spectral channels 350 and 680 a narrow- 
band interference signal. The fast- varying na- 
ture of the RFI signals added causes strong 
distortions of the reconstructed signal. For 
better visualization the plot shows a zoom-in. 
Fig. dZ] contains the quality indicators. 

channels and LOs containing RFI are more or 
less known (otherwise the fiagging would not 
have been possible). 

Implementing RFI flagging enables the 
analysis of the response of the LSFS to dif- 
ferent types of RFI. During our measure- 
ments we rarely encountered broad-band 
events, which last for only a second or less 
but affect several h undred spectral channels 
dWinkel et al.ll2007l ^. Fig.[20]and Fig.EUshow 
the result for affecting the 4th LO within spec- 
tral channels 200 to 400 — the reconstruction 
was successful when using our flagging scheme 
— it was not otherwise (not shown here) . The 
intensities of the broad-band signal are drawn 
from a power law but lie within 4 . . . 20o"rms- 
We added the interference onto each spec- 
trum of the 4th LO which would hardly be 
the case for a real observation (this type of 
RFI is rare). 
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Fig. 17. — The quality indicators (see Fig. [4] 
for the explanation of the symbols) are very 
sensitive to RFI signals. After subtracting a 
third-order baseline the signals noise level is 
increased by about 20%, while a is even in- 
creased by a factor of > 4. 

4. Computational efficiency 

When we started the analysis by imple- 
menting the LSFS within the C programming 
language we chose for the sake of simplicity 
the SVD algorithms delivered with the GNU 
Scientific Library (GSLj^]- They make use of 
the modified Golub-Reinsch algorithm. But 
in our case we can save a lot of computing 
time by using the fact that we have a sparse 
matrix. There exist a few libraries (mainly 
for FORTRAN) which use the Lanczos (SVD) 
algorithm for sparse matrices. We used the 
las2 routine from SVDPACKcH through the 
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Fig. 18. — Proper handling of RFI contami- 
nated data points allows reconstruction of the 
signal. In the signal domain residual RFI 
peaks remain, but have less amplitude and 
no measurable influence on their environment 
as this was seen in Fig.[16l We added three 
narrow-band interferences. The RFI signal 
in spectral channel 350 was added in all LO 
phases except for LO 1 and 2. The signal 
in spectral channel 500 affected every second 
LO phase while the signal in channel 680 was 
added in all phases. It turns out, that if an 
RFI signal is not persistent for a whole LO 
cycle the unaffected data points in the asso- 
ciated spectral bin can even be sufficient to 
reconstruct the signal without artifacts. The 
less LO phases are affected, the less impact of 
the RFI on the reconstructed signal is visible. 

interface library SVDLIBC|. Tabled] lists run- 
times of the pure SVD computation for differ- 
ent /, using both methods. The las2 algo- 
rithm is about an order of magnitude faster, 
which means it can calculate the SVD of a two 
times larger matrix within the same time (as 
the SVD computation scales as I^). There- 
fore, it is the preferred method for large val- 
ues of N ■ I. Note also, that the main mem- 
ory needed scales roughly as NI'^. For the 
largest of our problems (I = 2048, = 8) a 
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Fig. 19. — Using the flagging scheme to sup- 
press distortions by RFI signals provides noise 
level (RMS) values (bottom; see Fig. H] for the 
explanation of the symbols), which are only 
~ 14% higher and a value of a, which is ~ 25% 
higher than theoretically. 

1-GB-machine was barely sufficient using dou- 
ble precision arithmetic. 

As the SVD needs only to computed once 
per LO setup the more important contribu- 
tion to computing times needed is due to the 
LSFS calculation itself. Based on our experi- 
ence in many cases the convergence is reached 
after few (< 5) steps. When confronted with 
RFI etc. this increases up to 20 or more iter- 
ations until convergence. To account for this 
we monitor changes of the solution signal and 
break the iteration after the solution has sta- 
bilized. Based on different numbers of steps 
needed, the computation of the LSFS {N = 8, 
/ = 1024) takes ~ 0.05 s per iteration step on 
a modern desktop PC (2.0 GHz, x86). We al- 
ready used an optimized BLAS library which 
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Fig. 20. — LSFS for a broadband interference 
signal. As only one LO frequency is affected, 
the signal and bandpass could be well recov- 
ered. Fig. [21] shows the quality indicators 
which are only slightly increased compared to 
the undisturbed case. 
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Fig. 21. — Quality indicators (see Fig. U] for 
the explanation of the symbols) for spectra 
containing a broadband interference signal. 

makes use of SSE or equivalent features of 
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modern x86 cpu's meaning that there is prob- 
ably not much potential to speed up the com- 
putation of the LSFS. 

We also played around a little bit with the 
compiler extension OpenMF@ to parallelize 
the LSFS for use on multi-processor/core 
machines. This could improve the run-time 
by about 25% on a Dual-Xeon machine and 
about 15% on a Dual-Core processor. The 
maximum speed-up one could expect would 
be a factor of two. In fact, the LSFS computa- 
tion depends mainly on the multiplication of 
the (huge) correlation matrix with the input 
vector. Here the memory bandwidth has large 
impact on the overall speed which is possible. 



5. Summary 

In this paper we analyzed the statistical be- 
havior of the LSFS method as a function of in- 
tegration time as well as the robustness of this 
new method against various potential sources 
of errors as RFI signals and gain curve insta- 
bilities. It turned out that LSFS will provide 
very good solutions in most cases. However, in 
case of RFI the solution is strongly disturbed, 
rendering LSFS useless. We developed a flag- 
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Table 1: Computing times needed to calculate 
the SVD using different algorithms. 
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ging scheme which is able to deal with inter- 
ferences if there is a detection database con- 
taining accurate information where (in time 
and frequency) RFI signals were present. 

Mild bandpass instabilities are no prob- 
lem at all but very fast variations can cause 
moderate to severe distortions of the recon- 
structed gain curve. The latter, however, are 
far from realistic scenarios at modern radio 
telescopes, making LSFS the best choice even 
when confronted with ugly (but nearly time- 
independent) bandpass shapes. A strong ad- 
vantage of LSFS versus common frequency 
switching methods is that there is only a small 
frequency shift needed which results in much 
less bandpass variations at all. 

We also have shown that the LSFS will fail 
in presence of very strong emission lines as 
would be the case for example in galactic H I 
research due to the strong Milky Way H I line 
emission of the disk. Here the linear order ap- 
proximation is broken. We presented a possi- 
ble workaround by remapping the signal. 

While we could not work out a significant 
speed increase for the pure LSFS computa- 
tion, we at least propose the usage of a spe- 
cialized algorithm to compute the SVD of the 
correlation matrix. Such an algorithm turns 
out to operate an order of magnitude faster 
for sparse matrices, which in turn allows the 
computation of such a matrix for twice the 
number of spectral channels within the same 
computing time. 
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